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A non-perturbative Renormalization Group approach is used to calculate scaling functions for 
an 0(4) model in d — 3 dimensions in the presence of an external symmetry-breaking field. These 
scaling functions are important for the analysis of critical behavior in the 0(4) universality class. 
For example, the finite-temperature phase transition in QCD with two flavors is expected to fall 
into this class. Critical exponents are calculated in local potential approximation. Parameteriza- 
tions of the scaling functions for the order parameter and for the longitudinal susceptibility are 
given. Relations from universal scaling arguments between these scaling functions are investigated 
and confirmed. The expected asymptotic behavior of the scaling functions predicted by Griffiths 
is observed. Corrections to the scaling behavior at large values of the external field are studied 
qualitatively. These scaling corrections can become large, which might have implications for the 
scaling analysis of lattice QCD results. 
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I. INTRODUCTION 

The investigation of phase transitions and critical phenomena has been an important motivation in the development 
of the Renormalization Group (RG) [TJ [2J El HI El El El [SI [0] ? and it remains one of the core areas of the application 
of RG methods. Due to the universal behavior of different systems in the vicinity of a critical point, the results of 
relatively simple systems, such as spin models with an 0(N)-symmetry, play an important role in the analysis of phase 
transitions in much more complicated systems. 

Among the systems of interest is also QCD, where scaling behavior appears in various contexts. An example is the 
scaling behavior observed for the chiral phase boundary in the plane of temperature and number of massless quark 
flavors, near the critical number of flavors above which QCD is still asymptotic free but chirally symmetric in the 
infrared regime 112] . For two light quark flavors, the phase transition at finite temperature is expected to fall into 
the 0(4) universality class [TO]. While the dynamics of the non-Abelian gauge fields and quarks in QCD are difficult 
to describe in detail and require a non-pertubative treatment, universality arguments can be brought to bear on the 
analysis of the behavior near the phase transition. QCD lattice simulations remain indispensable as non-perturbative 
calculations in which all fcrmionic and gauge degrees of freedom are taken into account. However, since quark masses 
are large and explicitly break the symmetry, and since the simulation volumes are finite, the observation of critical 
behavior on the lattice is difficult. Strictly speaking, a phase transition cannot even occur in a finite volume. 

For this reason, the finite-volume behavior and the scaling behavior with an external symmetry-breaking field of a 
theory in the same universality class are of great interest for the analysis of the critical behavior in lattice simulations. 
Different methods can be used to study the scaling behavior in simpler model systems. However, since long-range 
fluctuations are essential for the critical behavior, even here a non-perturbative approach is needed. Renormalization 
Group methods are a natural choice for this problem. The RG approach used in this paper allows to study quantum 
field theories in an effective-potential framework in infinite as well as finite volume. It also allows the inclusion of an 
external symmetry-breaking field so that the dependence on the strength of such a field can be studied. A natural 
upper bound for the possible size of the field in such a calculation is set by the ultraviolet cutoff of the theory. For 
a spin model, such as the Ising model, such an external field corresponds to an external magnetic field H, whereas 
in QCD the bare quark masses fill the role of the symmetry-breaking field. In the case of a spontaneously broken 
symmetry, the external field controls the mass of the pseudo-Goldstone modes; in QCD, this is the mass of the pions. 
The renormalization group allows us to vary the external parameters of a given theory, such as volume size or external 
symmetry breaking, over a wide range in parameter space, from a rather deformed theory to the theory of physical 
interest. This makes it possible to close a gap between scaling behavior in a finite volume and with a given external 
field and the scaling behavior in the limit of interest. 

Even though current lattice calculations are performed in large volumes and close to the physical value of the pion 
mass, it remains difficult to extract the order of the phase transition, and this question is still not conclusively settled 
[131 ITU fT5] . Usually, only the known critical exponents are used as input for a finite-size scaling analysis of lattice 
data. Including information about the scaling function as well would increase the power of the analysis HH [J71 [HI [19] . 
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The question of the size of the scaling region also remains open, at least for staggered fermions current volumes might 
still be to small to observe finite-size scaling [20] . 

The main results of the present paper are accurate determinations of the scaling functions of an 0(4) symmetric 
model in three dimensions over a wide parameter range for the external symmetry breaking field. Scaling functions for 
the order parameter and the longitudinal susceptibility arc calculated, and the relations between these functions are 
investigated. We also investigate qualitatively scaling corrections for large values of the external field. The discussion 
of finite-size scaling functions from a finite-volume calculation |21l I22j is the subject of a forthcoming paper. The 
discussion of other 0(N)-symmetric models will also be postponed to a future publication. 

The paper is organized as follows: In Sec. [TTJ, we provide an overview over the scaling behavior in the three- 
dimensional 0(4) model, known results for the scaling function and relations between the different forms of the 
scaling function. The RG method is introduced in Sec. |III| and technical details and the setup of our calculation 
are discussed. We determine the critical exponents for our calculation in Sec. |IV| to ensure a consistent analysis 
of our results. A comparison to values from other methods allows to assess the systematic errors inherent in our 
calculation. Our main results for the scaling functions are presented in Sec.[V]for the order parameter, and Sec.|VI|for 
the susceptibility. Sec. |VII| contains a brief discussion of the masses of the longitudinal and transverse fluctuations, 



m a an m v , and we present our concluding remarks in Sec. VIII 



II. SCALING IN THE 0(4) MODEL IN d = 3 

Scaling behavior can be observed in the vicinity of a critical point. Critical behavior is characterized by a diverging 
correlation length £ — * oo. Close to a critical point, there is no relevant length scale due to the critical long- 
range fluctuations, and such a system is invariant under a change of length scale. As a consequence, the behavior 
of thermodynamic observables in the vicinity of the critical point is characterized by critical exponents which are 
universal for systems in the same universality class. 

The critical behavior is governed by the singular part of the free energy, and the behavior of thermodynamic 
variables can be derived from the free energy. We consider an 0(4) model with an external symmetry-breaking field 
H and the temperature T as the two relevant coupling constants. The invariance of the singular free energy density 
under a rescaling of the length with a factor I close to the critical point can be expressed as 

f s (t,h) = r d f s (&%^h). (i) 

Possible corrections to the scaling behavior due to irrelevant couplings are neglected. The dimensionless couplings t 
and h are defined as 

r£ rp TJ 

t=— - and h=—, (2) 

Jo -Ho 

and the critical point is located at (T, H) = (T C) 0) or (t, h) = (0, 0). The critical exponents can be expressed in terms 
of y t and yh and are connected by a number of scaling laws, so that there are only two independent ones: 

1 v 

yt = -, yh 



V 



/3<T 



7 =(2- 77 ) 2 y, 7 = /3(«5-l) ! /3 = l -{d - 2 + rj)v, ud=0(l + S). (3) 

The critical exponent v governs the behavior of the correlation length £ and the critical exponent 7 governs the 
behavior of the (longitudinal) susceptibility \ as the critical temperature is approached, 

£~i*r, x~i*r- (4) 

The critical exponents (3 and S describe the behavior of the order parameter M as the critical point is approached 
from either the h- or the ^-direction: 

M(t = 0,h) = h 1/s , M(t, h = 0) = {-t) for t < 0. (5) 

In the usual convention, which we will use in our analysis, the normalization constants Hq and To are determined by 
these two relations. 
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Starting from the scaling form of the free energy Q, the scaling behavior of the order parameter M and the 
longitudinal susceptibility \ can be derived using their thermodynamic definitions 

M = -§ (6) 
dM 

x = dW (7) 

by choosing the scaling factor I in an appropriate way. The most intuitive and most commonly used form of the 
scaling functions is obtained by choosing £ Vh h = 1. The free energy density then becomes a function of h and the 
scaling variable z — t/(h 1 ^^) only, and the order parameter can be expressed as 

M = h^f(z)^^L = f(z), z^j^y (8) 

where f(z) is a universal scaling function. With the normalization ([5|, f(z) has the properties /(0) = 1 and f(z) — > 
{—z)p for z — > — oo. 

From the definition of the longitudinal susceptibility x, it follows that its behavior can also be described in terms 
of a scaling function f x (z), which in turn is determined from the scaling function of the order parameter and its 
derivative: 



X H 5 



&H h 1 - 1 / 5 X = f x (z) = i 



(9) 



The scaling function for the susceptibility is therefore completely determined if one knows the critical exponents and 
the scaling function for the order parameter. This relation presents an additional, very non-trivial test for the scaling 
behavior. 

While the scaling function f(z) expressed in terms of the scaling variable z is convenient for visualizing the behavior 
of the order parameter, it is not the most stringent test of the scaling behavior as it is relatively insensitive to scaling 
violations. In addition, we find that it is not very convenient for parameterizing the scaling behavior. Better for this 
purpose is the Widom-Griffiths parameterization j23j EI] of the equation of state as y = y(x) in terms of the scaling 
variables 

t h . 

y^Tis- (io) 



MV/3 ' y M 5 ' 

(For definiteness, we will take H > and M > throughout.) This parameterization can be obtained from the scaling 
form of the free energy (JT|) by taking £ y *t = 1. It is equivalent to the first one, and the two parameterizations are 
related by 

With the normalization ([5]), the function y = y(x) satisfies the normalization conditions y(0) = 1 and y{— 1) = 0. 

Assuming that the equation of state in the form H = H(M, T) is analytic everywhere away from the first-order line 
H = 0, T < T c and the critical point H = 0, T = T c , Griffiths showed [21] that for H > 0, T > T c and M > M a > 
for some value Mq, i.e. x > 0, the equation of state can be expressed in a convergent expansion as 

oo 

y{ X ) = ^ CbI 1-*- 1 ) =^( Cl + C 2 X~ 2/? + C3X~ 4/3 +...). (12) 

n=l 

For asymptotically large values of the scaling variable x, this expression allows us to determine the critical exponent 
7 from the leading term and thus provides an additional check for the consistency of the scaling behavior. It also 
suggests a parameterization of the equation of state y = y{x) for large x. 

As for the first parameterization, the scaling behavior of the susceptibility can be expressed in terms of the scaling 
function y = y(x), its derivative and the critical exponents. One finds 

X = qjj = ms^ij^ ( Sy(x) - pxy'(x) ) 4* [HoM'-'x]- 1 =Sy(x)--xy'{x). (13) 

We consider later the inverse rescaled expression [H () M s ~ 1 x\ 1 to avoid the divergence of the universal result which 
appears at x — in this parameterization. 
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With Griffiths' expression (12), and using the scaling law 7 = (3(5 — 1), one finds for the rescaled susceptibility for 
large x the expression 

[HqM^x] = S (y(x) - j- § x y'(xf) = x~< ( Cl + 3c 2 x^ + 5csX~^ + . . .) 

OO 

= ^c n (2n-l) a; T- 2 «"- 1 ), (14) 

n=l 



with the known coefficients c n from (|12j) . Remarkably, this means that for large x the leading-order terms from the 
scaling function of the order parameter y(x) — Cix 1 + . . . and the scaling function of the inverse rescaled susceptibility, 
Sy(x) — jjXy'(x) — cix 1 + . . . should be exactly the same, provided the scaling laws for the critical exponents hold. 
The transverse susceptibility 

XT = ~ (15) 

also satisfies a universal scaling relation. With eq. ([8]), one finds that the scaling function is given by the one for the 
order parameter, after \t nas been rescaled in the same way as the longitudinal susceptibility: 



h 1 ' 5 - 1 (^) = ^ W /W * H h^ XT = f(z). (16) 



1 , 1 f M 



Using the Widom-Griffiths scaling form, one finds likewise that the rescaled susceptibility is given by the scaling 
function y(x): 

XT = J^l] * [HoM^xtV 1 = y. (17) 

Since we only recover the scaling function of the order parameter from these relations, the transverse susceptibility 
does not provide any additional test of the scaling behavior. In our implementation of the model, the relation (151 
follows from the determination of the order parameter from the minimum of the potential and is explicitly satisfied 
by construction. For this reason we do not analyze the scaling relation for \t separately. 

There are few results for the scaling functions in the literature. The equation of state in the form y = y(x) was 
calculated in the e-expansion to 0(e 2 ) by Brezin, Wallace and Wilson [55], and a parameterization of the result close 
to the coexistence point x = — 1 (corresponding to the first-order line H = 0,T < T c ) is provided by Wallace and Zia 
|26J. Monte-Carlo simulations of the 0(4)-symmetric spin model on the lattice were used by Toussaint [16 to obtain 
the scaling function f(z), and by Engels and Mendes (27] to obtain a parameterization for the equation of state y(x) 
(see also |18j). and by Cucchieri and Mendes for a comparison in a parametric form |28j . 

The result obtained from the e-expansion |25j for N = 4 in d = 3 becomes after expansion around x = —1 |26j 

y(x) = (l + a;){l + e[Alog(l + x) + B] + e 2 [C , log 2 (l + a;) + Z)log(l + a;)+£;]+C'(e 3 )} (18) 
with A = 0.125, B — —0.0270494, C = -0.015625, D = 0.148902, E = 0.138468. 

This expression is valid for (1 + x) <^ 1. In [26J it is argued that the expansion around x = — 1 can be inverted to 
give an equation of state x — x(y) in the form 

x{y) + l = hy + 5 2 y W2 + diy 2 + d 2 y 2 ~ e/2 + d 3 y 2 ~ e + . . . , (19) 

where the coefficients are determined only to 0(e) while the original expansion is known to C(e 2 ). For N = 4 and 
d = 3 (e = 1) the perturbative result from [26] is 

x(y) + 1 = c 2 y 1/2 + (ci + d 3 )y with (20) 

5 2 = 0.250 + 0.280e = 0.530, c x = 0.750 - 0.253e = 0.497, d 3 = 0.03125, c x + d 3 = 0.528. 

The normalization condition y(0) — 1 is only satisfied to about 10% (y(0) = 0.92 for the inverted expansion), and 
the expansion is very accurate only close to the expansion point x = — 1. In general, it is known that the e-expansion 
satisfies the large- a; behavior from Griffiths' expansion (12 1 only order by order in e, but not explicitely [25], and thus 
even the full expression is not expected to be accurate for i» 1. 
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While critical exponents from the e-expansion are available to higher order in e [55], for consistency we use in the 
comparisons involving the scaling function the critical exponents calculated to 0(e 2 ) for the 0(4) model in d = 3 [25] : 

1 f\ 1 QQ 1 H7 

/9= — sa 0.3902, 7 = « 1.3854, 6= « 4.458. 

41 '96 24 

These values satisfy the scaling law ") = /3(<5 — 1) to within 3% 

We observe that an ansatz y(x) = c (l + x) T reproduces the first terms of the e-expansion around x = —1, if we 
assume that the exponent can be expanded as r = 1 + ne + . . . in e: 

y{x) = = cq(1 + x) 1+Tl£ +- = cq(1 + x)(l + nelog(l + a:) + . . .)• (21) 

For this reason we use such an ansatz with the exponent r as a free parameter to fit our results, and we find that 
such a parameterization works surprisingly well. In addition, for large x S> 1 it coincides with the leading term of 
Griffiths' expression if r = 7. 

From Monte-Carlo lattice simulations of an 0(4) spin model, the parameterization for the equation of state from 
Engcls and Mendes |27j is also provided in the inverted form x — x{y). For small values of x, the parameterization 
x s (y) is taken from the e-expansion (21 1, with coefficients determined non-perturbatively from a fit to the lattice 
results PHGI7]: 

x s (y) + 1 = c 2 y 1/2 + (ci + d 3 )y + d 2 y 3/2 with (22) 

~c 2 = 0.674(08), ci + d 3 = 0.345(12), d 2 = -0.023(5). 



For large values of x, the first two terms of an inverted form of Griffiths' expression ( 12 1 are used [18] : 

xi(y) = ay 1 /'' + by^~ 2 ^h with 
a = 1.084(6), 6= -0.994(109). (23) 

The values of the critical exponents are here 

(3 = 0.380, (5 = 4.86 7= 1.4668, v = 0.7423. 

The authors also formulate an interpolating expression for the equation of state 

3 3 

x{y) = x s (y)^^+x l (y)^-^, y = 10.0 (24) 

which describes the results well for a range of at least — 1 < x < 30 and < y < 150. Inspired by this idea, we will 
use a similar interpolation for the equation of state y = y(x) to describe our results. 

In |28j . in addition to a parametric description, the results for the scaling function is also given in the above form. 
The coefficients found there are 

5 2 = 0.746(3), ci + d 3 =0.19(1), d 2 = 0.061(8) (25) 
for the small- a; region (x < 1.5), and 

a =1.07(1), 5= -0.95(3) (26) 

for the large- x region. The same expression is used for an interpolation. The differences to the earlier results are 
attributed to slightly different values of the critical exponents used in the calculation. 



III. RENORMALIZATION GROUP APPROACH TO SCALING 



In this section we discuss our Renormalization Group (RG) approach to scaling behavior in the presence of an 
external symmetry-breaking field and the derivation of the proper-time flow equations for the 0(4)-potential in 
infinite volume. We also give a detailed discussion of the approximations that we have used throughout this work. 
Reviews of and introductions to functional RG methods can be found in e.g. [30] |3TJ [32] [33] [34] [35l [36] 157] . The 
effective action of the 0(4) model in d space-time dimension is given by 

T[cl>] = J d d x {\z^dj) 2 + . (27) 
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Here we have neglected possible kinetic terms of higher order such as Y^^d^) 2 . The components of the vector <f> 
are labeled according the role that the corresponding fields are playing in the spontaneously broken regime, cp- = 
(a, 7T 1 , 7r 2 , 7r 2 ). We choose the first component to be the radial mode in the regime where the ground state of the 
theory is not symmetric under 0(4) transformations: 

(<?) =0o T = (to, 0,0,0) (28) 

The potential U((f)) in its present form depends only on <f> 2 = (t 2 + tt 2 . As discussed in [38, 39, 40 such an ansatz for the 
potential is not appropriate if one is interested in a study of phase transitions in a finite volume, since strictly speaking 
spontaneous symmetry breaking does not occur in a finite volume. Therefore the presence of an non- vanishing external 
source in the ansatz for the effective action is indispensable for our study. Moreover we do not solve the RG flow for 
the full potential but expand the potential in a Taylor series around the vacuum expectation value and solve then the 
RG flow for the expansion coefficients of this series. Before we discuss our ansatz for the potential in more detail, we 
briefly sketch the derivation of proper-time RG equations in general. 

In a Gaussian approximation, we can perform the functional integration of the bosonic fields and obtain the one-loop 
effective action for the scalar fields <fi, 

m = T Auv [0] + ^Tr log (vf [0]) (29) 

where T B [(f)] is the inverse two-point function evaluated at the background-field and T\ contains the initial values 
of the RG flow at ultraviolet (UV) scale A. From now on, we neglect a possible space dependence of the expectation 
value and take the wave-function renormalization to be constant, Z$ = 1. Therefore the anomalous dimension of 
the scalar field is zero, r\ — 0. We have to keep this approximation in mind for the interpretation of our results, since 
the anomalous dimension of the scalar fields influences the critical exponents which in turn are the key ingredients for 
finite-size scaling. However, the scalar anomalous dimension is small compared to one, see e. g. [41] . and therefore 
our approximation to neglect the running of the wave-function renormalization is well-justified. 



In order to regulate the infrared (IR) divergences in Eq. ( 29 ) we use the Schwingcr proper-time representation of 
the logarithm and introduce an IR cutoff function 1 / a (rfc 2 ), where the variable r denotes Schwinger's proper time 
and k is a cutoff scale which has mass dimension one. The derivative of the cutoff function with respect to the scale 
k is given by (see e. g. Refs. [71 14"3"]) 

In the following we choose a — ^. The relation between the so-called proper-time RG and the Functional RG has 
been worked out in detail in [JH 351 H5J. As was found in Refs. [UJ 28], the flow equation which results from this 
choice for the proper-time cutoff function coincides with the flow equation obtained in the Functional RG framework 
with an optimized cutoff. 



/ 9\ . . 

The inverse two-point function [</>] in Eq. (29 1 depends on the second derivatives of the effective potential U, 



The eigenvalues of this matrix evaluated at the minimum of the potential are the masses of the fields. By replacing 
the bare masses and couplings in the inverse two-point functions with the scale-dependent quantities, we obtain the 
so-called renormalization group improved flow equation for the effective potential in infinite volume (see also e. g. 
01171 HE US]): 

d - (fc 2 ) d/2+1 l/3 1 \ 

k dk Uk ^ = (4vr)<V 2 r(d/2+l) \k 2 + Ml k (<P) + k 2 + M 2 a k {<P) J ' (32) 

The quantities M a and M n are the eigenvalues of the the second-derivative matrix of the potential. Note that these 
qunatities still depend on the background field (/). 



On the one hand physical quantities calculated from the RG flow should not depend on the choice of the cutoff function in the limit 
k — > 0. On the other hand, a study of the regulator dependence of the results in this limit allows us to check the quality of our 
truncations. A detailed study of the dependence of our results on the choice of the regulator function is in progress 42 . 
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For a study of critical exponents and scaling it is convenient to deal with dimensionless quantities rather than 
dimensionful quantities. Therefore we introduce the dimensionless potential u, the dimensionless masses m a and 
as well as the dimensionless field ip by 



u k (<p) = k d U k {4>) m? fc fo>) = AT 2 M? fc (0) and ip = k~ 



Applying these definitions to the flow equation ( 32 1 , we obtain 

1 



d t ut{ip) 



-du(tp) 



1 



( 4?r )d/2 r( d / 2 + 1) \\ + ml{ip) 1 + ml{<p) 



(33) 



(34) 



where the dimensionless flow variable t is given by t = ln(fc/A). By integrating the flow equation (either the dimen- 
sionless or the dimensionful formulation) from the UV scale A to k — *■ 0, we obtain an effective potential in which 
quantum corrections from all scales have been systematically included. 

We now discuss the ansatz for the effective potential £/&(</>). Since we are interested in a study of phase transitions 
with an external symmetry-breaking field, we have to introduce a corresponding linear term in the field into our ansatz 
for the effective action. In order to solve the RG flow for the effective potential U, we expand the potential in local 
n-point couplings around its minimum ao(k) 



U k (a, n 2 ) = a (k) + ai (k){a 2 + tt 2 - a {k) 2 ) + a 2 {k){a 2 + tt 2 - a {k) 2 ) 2 + . 



Ha 



(35) 



where H is the fixed, external symmetry-breaking field. The minimum ao(k) is the order parameter of the system 
which can be identified with the magnetization M in an Ising model, 



M = a . 



(36) 



Since we have absorbed the symmetry-breaking linear term in the effective action into the ansatz for the potential, 
Uk depends now on the fields a and 7r separately. The condition 



dU k (a,Tr 2 ) 



(37) 



ensures that we are expanding around the actual physical minimum. From Eq. ( 37 ) , we find that the RG flow of the 
coupling a\ and the minimum a are related by the condition 



2ai(k)a (k) = H. 



(38) 



This condition keeps the minimum at (cr,w) = (<7o(fc),0). The flow equation of the minimum is thus related to the 
flow of the coupling a± in a simple way. Taking the derivate with respect to k of Eq. ( 35 1 , we obtain 



dU k (a,7T 2 ) _ u dao(k) 
k ~ ~dk 



2ai(k)ao(k)k 



da Q (k) 
dk 



ok 



(39) 



Note that the derivative of the potential with respect to the regulator scale does not contain any contributions linear in 
the fields. Of particular interest for our studies in this work is, apart from the order parameter qp (fc), the longitudinal 
as well as the transversal susceptibility. Taking the total derivative of the minimum condition ( |37[ ) with respect to 
H, one finds that the longitudinal susceptibility is inverse proportional to inverse mass squared of the radial mode: 



1 



1 



XL 



2a\ + 4a 2 a a 2 



(40) 



The transverse susceptibility is defined as the ratio of the order parameter M = a (magnetization) to the external 
field H and is inverse proportional to the mass of the (pseudo-) Goldstone particles 



XT 



M 



1 

2ai 



(41) 



The RG flow equation for the couplings Oj can now be obtained straightforwardly by expanding the flow equation 
( |32[ ) around the minimum ao(k) and projecting it on Eq. (39). The RG flow of the minimum of the potential is 
determined by Eq. ( 38 1 . Apparently, this procedure for a derivation of the flow of the potential u results in an infinite 
set of flow equations for the couplings a*. In order to solve the set of equations for the couplings, we have to truncate 
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our ansatz (35 1 for the potential. In the following, we include fluctuations around the minimum up to eighth order 
in the fields, i. e. we keep track of the running of the couplings a±, a 2l 03, and 04. The resulting finite set of 
coupled first-order differential equations are then solved numerically. We have checked that our results do not change 
significantly by including couplings of higher order. We have also checked numerically that our results coincide with 
earlier results (e. g. [HUH]) in the limit of vanishing external source H. 

For our calculations, we have chosen a cutoff scale A = 1000 MeV, which is comparable to a typical lattice cutoff 
(i/a « 1500 MeV with a being the lattice spacing) in current finite-temperature lattice simulations. The initial 
values for the differential equations of the couplings are given by the effective potential at the UV scale in the form 
U\(a, 7? 2 ) = a 2 (A)(0 2 — 0q(A)) 2 — Ha, where we choose 02 (A) /A = 0.025. The flow of a\ is determined by the flow of 
the minimum <7o(fc). In e? = 3 dimensions, the initial value of the minimum of the symmetric potential at the cutoff 
scale 0o(A) serves as proxy for the temperature, and we assume that an expansion (tpo(A) — </>Q rltlcal (A)) ~ (T — T c ) is 
possible. This is to be seen in analogy to lattice simulations of spin models, where the temperature is absorbed into 
the spin-spin-coupling, which can then be tuned such that the system is at its critical point. Similarly, in Landau 
effective theory it is assumed that the couplings have an expansion in T — T c around the critical values. 



IV. CRITICAL EXPONENTS 



Critical exponents for the 0(4)-model in d = 3 have been determined by many different methods to a high degree 
of accuracy. Theoretical results are available from perturbative field-theoretical calculations 29, 501 [FT], from lattice 
Monte-Carlo simulations of spin models [521 155] and 4 theory [53], and from non-perturbative RG calculations 
[4"T1 |4"T1 l4"51 l4"9l 1551 [55] , (See e.g. [57] for an overview.) It is not the purpose of the present paper to add to this list, 
in particular since our results coincide with [48l and since we work in the local-potential approximation where the 
anomalous dimension vanishes, 77 — 0. More accurate results from functional RG methods beyond this approximation 
are available [41, 49, 55]. But in order to present a consistent evaluation and to make contact with these calculations, 
we determine critical exponents from our results and use these values in the analysis. 

Since we are working in d = 3, temperature cannot be defined field-theoretically. However, the relevant coupling 
for approaching the critical point is the initial value of the expectation value of the field <fi a t the ultraviolet scale A, 
</>o (A), and we assume that an expansion (</>o(A) — </>Q rltlcal (A)) ~ (T — T c ) exists. 

With our choice of scale A = 1000 MeV and initial parameters, we determine the non-universal critical temperature 
of our model in the absence of an external symmetry-breaking field to be 

T c /A 1/2 = 0g ritical (A) = 13.682 368 165 072 75(1) MeV 1/2 . 

(We find that this accuracy is actually necessary in the choice of the initial conditions to see the fixed point behavior 
clearly in the RG flow.) 

We determine the critical exponents (3 and 5 from the order parameter M, and the critical exponent v from the 
correlation length £ — l/m a in the phase with restored O(N) symmetry. The critical critical exponent (3 can be 
determined from the order parameter at H = from 

logM = /91og(-t) = -/31og(T ) + /?log(T c - T), (42) 

by regarding log M as a function of the variable log(T — T c ) and either fitting this directly as a linear expression 
in a region close to T c , or by numerically taking the derivative of the function in the limit T — > T c . From an RG 
point of view, the critical exponents arc defined from the eigenvalues of the linearized flow equations at the critical 
point. Methodically, it is therefore more sound to use the limit of the derivative for the determination. We do indeed 
find good agreement with the results of [48j 158] . where the critical exponents were determined in local potential 
approximation from the same RG flow equations by diagonalizing the stability matrix at the fixed point. 

In order to estimate how well the critical exponents can be determined from our calculation, we also fit the linear 
expression below T c . We repeat the determination from the derivative and from the fit in the presence of a small 
external field H — 1.0 x 10~ 9 MeV 5 / 2 . The results are given in Tab. [i] While the original determination from the 
derivative gives the most exact value, we use the results from the other determinations in order to estimate an error 
for this determination. We observe that the values from a fit tend to be larger for (3 and smaller for u, compared to 
the values from the derivative. The values for the critical exponents (3 and v satisfy the scaling law (3 — \{d — 2 + r\)v 
with 77 = and d — 3 with better than 0.2%. 

The critical exponent 5 can be determined from the order parameter M at the critical temperature T c , taken as a 
function of H 

log M = - 1 log H + i log H (43) 
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H = 





H = 1.0 x 10 


" 9 MeV 5/2 




derivative 


fit 


derivative 


fit 


quoted 


P 0.4030 


0.4051 


0.4033 


0.4108 


0.4030(30) 


v 0.8053 


0.7953 


0.8007 


0.7920 


0.8053(60) 



TABLE I: Critical exponents /3 and v determined from our results. For the evaluation, we use the quoted values. The errors 
reflect the uncertainty of the determination, in addition there is a systematic error for the calculation which is not quoted and 
which cannot be determined from the present calculation taken by itself. 



derivative fit vd/ft—l (d + 2 — fl)/(d — 2 + rj) quoted 
5 4.9727 4.8409 4.995(90) 5.0000 4.973(30) 



TABLE II: The critical exponent 8, determined directly from the derivative and from a fit to the order parameter at the 
critical temperature as a function of H, and indirectly from the other critical exponents by means of the scaling laws. In our 
calculation, the anomalous dimension r\ — 0, and thus S = 5 is theoretically expected. 



again by either calculating the derivative in the limit H — > or by fitting this as a linear expression in logH. We 
expect the result from the limit of the derivative to be more reliable, but also calculate S from a fit to estimate the 
accuracy of our determination. The results are given in Tab. [TTJ In contrast to our direct determination, in most 
functional RG calculations |4T[ 149] . the critical exponent S is calculated from the anomalous dimension rj with the 
scaling law 

a — 2 + rj 

or determined from the asymptotic behavior of the effective potential. Since rj = in the local potential approximation, 
we expect 6 = 5 from the scaling law. In addition, we can also use the values for the critical exponents (3 and v to 
calculate 5 with the scaling law dv — /3(1 + ft). The results are also given in the table and are used in the estimation 
of the error of the determination. 

The total error for the values of the critical exponents is larger than the uncertainty from the determination which 
we give explicitly. It is dominated by the larger systematic error from the approximations that are necessary to solve 
the RG flow equations. This systematic error cannot be estimated by looking at the current calculation in isolation. It 
is generally difficult to assess these systematic errors in non-perturbative functional RG calculation, but the errors due 
to the necessary truncations can be estimated by comparing different RG cutoff schemes for the same approximation, 
and by systematically improving the approximations. The two main approximations in the present calculation are a 
restriction to only a finite number of n-point couplings in the effective potential, and a truncation in the momentum 
dependence to local couplings. 

The restriction to a finite number of n-point couplings is not very severe, and its effects can be estimated rather 
reliably by including a larger number of couplings. The results in [HI 08] show that the convergence is relatively fast, 
and the inclusion of a small number of n-point couplings is usually sufficient in the case of O(N) models. Our own 
observations confirm these findings. 

The effects of the truncation in the momentum dependence of the couplings (derivative expansion) are more difficult 
to assess. In the current calculation, the vanishing anomalous dimension rj, which is due to the restriction to local 
couplings, affects all coefficients, as can be seen from the scaling laws. For this reason we do not expect to obtain 
values for the critical exponents comparable to those from simulations, where the anomalous dimension is non-zero. 

For our results, the effects from this truncation can be best estimated by comparing to other results from functional 



RG calculations. We list a number of results in Tab. Ill We expect our results to coincide with those of [48], since this 
determination of the critical exponents from the fixed point is performed in the same RG cutoff scheme and in local 
potential approximation, and we do indeed observe good agreement. The FRG calculations [UJUiOES] au include a 
wave function renormalization, which leads to a non-zero anomalous dimension. As one can see from the table, the 
resulting critical exponents are much closer to the values observed in lattice Monte-Carlo simulations. The differences 
to our results are a measure for the systematic error due to the local potential approximation. 

Nevertheless, the values obtained here are the ones needed for a scaling analysis of the results from our calculation, 
since the effects that are responsible for the non-zero anomalous dimension are absent in the effective potential and 
our observables as well. For this reason we must use these values for a consistent evaluation of the scaling behavior. 

Once the values for the critical exponents (3 and 8 are known, we use these values to determine the normalization 
constants T and H from fits to the results. For To, we do this by fitting with fixed (3 to the order parameter for 
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method 


;/ 


P 


V 


S 


R. Guida, and J. Zinn- Justin 


[29] FT 


0.741(6) 


0.3830(45) 


0.0350(45) 


4.797(25)* 


R. Guida, and J. Zinn- Justin 


[29] e-exp 


0.737(8) 


0.3820(25) 


0.036(4) 


4.792(22)* 


K. Kanaya, and S. Kaya 


[52] MC 


0.7479(90) 0.3836(46) 


0.0254(38)* 4.851(22) 


H. G. Ballesteros et al. 


[S3] MC 


0.7525(10) 0.3907(10)* 0.0384(12) 


4.778(8)* 


M. Hasenbusch 


[SI] MC <f> 4 


0.749(2) 


0.388(2)* 


0.0365(10) 


4.789(6)* 


G. v. Gersdorff, and C. Wetterich gS] FRG 


0.739 


0.387* 


0.047 


4.73* 


N. Tetradis, and C. Wetterich 


[41] FRG 


0.791 


0.409 


0.034 


4.80* 


O. Bohr et al. 


g9] FRG 


0.78 


0.40 


0.037 


4.80 


D. F. Litim, and J. M. Pawlowski 


g8] FRG 


0.8043 


0.4022* 




5.00* 


our result 


FRG 


0.8053(60) 0.4030(30) 




4.973(30) 



TABLE III: Comparison of critical exponents for the 0(4) model in d — 3 from different methods. Values with an asterisk are 
not obtained independently, but calculated from other exponents with the scaling laws. This table is by no means comprehensive, 
see e.g. [S7J for a more complete listing. In addition to the functional RG (FRG) results, we list results from lattice Monte-Carlo 
(MC) calculations, from perturbative field-theoretical (FT) calculations in d — 3, and from the e-expansion. 



small values of the external symmetry breaking field H — 1.0 x 10~ 10 MeV 5 / 2 and H = 1.0 x 10~ n MeV 5 / 2 . For H , 
we fit results obtained exactly at T Cl using the fixed value of S. We find for the non-universal normalization constants 
in our calculation 

T /A 1/2 = 0.014916(5) MeV 1/2 , H a = 6.032(10) MeV 5/2 . (45) 
In the following, we use these values in our analysis. 



V. ORDER PARAMETER 



For the determination of the scaling function, we have calculated results for the order parameter M, the longitudinal 
susceptibility (or mass m a of the radial excitations) , and the transverse susceptibility (or mass of the transverse 
excitations), over a wide range of values for T and H. For H, the values span a range of seven orders of magnitude, 
from H — 1.0 x 10~ 4 MeV 5 / 2 to H = 1.0 x 10 3 MeV 5 / 2 . For each value, we have chosen appropriate ranges for T 
around T c , so that we always cover a similar range of values for the scaling variable z = —10 ... 10, or z = —30 ... 30 
for the larger values of H . 

In this way we can determine the scaling functions with high accuracy from the results with small symmetry- 
breaking fields, where scaling corrections remain small. At the same time, we can assess scaling violations from the 
results at large fields and determine where these corrections become large. This is relevant for the application of a 
scaling analysis to a physical system, and our range of values extends beyond the currently available results from spin 
model lattice simulations. 

In the following, we use the normalization constants and critical exponents determined in section IV Since we work 
in an approximation with local couplings the anomalous dimension in our calculation is zero, -q — 0. This introduces 
a systematic error into critical exponents and observables which cannot be estimated from our calculation taken by 
itself. Within our approach, the systematic error can only be estimated by comparing to other RG calculations which 
go beyond this approximation or which use a different cutoff scheme. Work on a comparison of different RG schemes 
and on an improvement of the approximation are in progress. It should be kept in mind that for now this systematic 
error remains unquantified in our results. Of course, it can also be estimated by comparing directly to spin model 
lattice simulations. 

However, despite the observed deviation from the lattice in the current approximation, we will see below that 
the results are internally consistent: The scaling laws among the critical exponents are satisfied very well, and the 
expected scaling behavior is observed to an impressive accuracy. The same values of the critical exponents are also 
recovered from the asymptotic behavior of the scaling functions, as expected. In addition, the scaling functions for 
the order parameter and for the susceptibility are well consistent with each other, as we will demonstrate below. 
These observations serve as a cross-check of our results and confirm that the critical exponents have been determined 
correctly for our calculation. 

For our determination of the scaling function, we only consider the leading-order scaling behavior and do not 
quantify the observed corrections to scaling. For the purpose of judging the quality of a fit to the leading-order 
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s-range 


c r 


x 2 


# d.o.f. 


X 2 /# d.o.f. 


-1< x < 1 


0.9928(14) 1.6712(27) 


568 


838 


0.6778 


-1< x < 5.3 x 10 4 


0.9928(21) 1.5941(8) 


3682 


1551 


2.374 



TABLE IV: Fit to the rescaled order parameter for H = 1.0xl0~ 4 , l.OxlO -3 , 1.0xl0~ 2 MeV 5/2 with the ansatz y(x) = c(l+x) T . 
The scaling deviation is estimated from the spread of the rescaled susceptibility. Shown are separate fits for the range — 1 < x < 1 
and for the global fit — 1 < x < 5.3 x 10 4 . 




FIG. 1: Scaled RG results y(x) for H = 1.0 x 10~ 4 , 1.0 x 10~ 3 , and 1.0 x 10" 2 MeV 5/2 . Errors are estimated from the spread 
of the data points after rescaling and are due only to scaling corrections. Not all data points are displayed. Shown are further 
fits to the small- a; values x < 1 (solid black line, small- a; fit) and to all x-values up to approximately x — 5.3 x 10 4 (dashed blue 
line, asymptotic fit). 

scaling behavior, we treat scaling deviations due to a small variation of the external field as an error. This is suitable 
for weighing how well a single scaling function can fit the results after rescaling. We estimate the scaling corrections 
by comparing the rescaled results for different values of H. 

In the following, we first determine a fit for the scaling function in Widom-Grimths parameterization. We found 
this to be the most convenient parametrization, and also the one most sensitive to deviations from scaling, so that it 
provides the most stringent test. We will then compare the results to those from the lattice and the e-expansion. After 
translating the results into the more intuitive scaling form f(z) as a function of the scaling variable z, we compare 
our parameterization to results at large values of the external field and demonstrate deviations from scaling. 

A. Scaling function in Widom-GrifHths parameterization 

1. Phenomenological ansatz for the scaling function 

We start our analysis with the results for H = 1.0 X 1CT 4 MeV 5 / 2 , H = 1.0 X 10" 3 MeV 5 / 2 , and H = 1.0 X 10~ 2 
MeV 5 / 2 , since corrections to the leading-order scaling behavior are small in this range. This allows us to obtain a 
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2.0x10 



1.5x10 



y l.oxio 7 - 



— asymptotic fit 
° H = 1.0 10" 4 MeV 5 ' 2 
#=1.0 10" 3 MeV 5 ' 2 
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FIG. 2: Scaled RG results y{x) for H = 1.0 x 10~ 4 , 1.0 x 10 -3 , and 1.0 x 10~ 2 MeV 5/2 . Errors are estimated from the spread 
of the data points after rescaling (so they are mainly due to the appearance of scaling corrections). Not all data points are 
displayed. The fit to all a;- values up to approximately x — 5.3 x 10 4 is also shown (black dashed line, asymptotic fit). 



simple but quite decent parameterization of the scaling function y — y{x). For the fitting procedure, we estimate 
errors due to corrections to scaling by comparing the results after rescaling for the different values of H. 

For small values of x close to the coexistence point x = —1, it is natural to expand the equation of state in terms 



of (1 + x). Inspired by eq. (21), we fit the phenomenological expression 

y(x) = c(l + x) T (46) 
to our results, which provides a surprisingly good fit with only two parameters. Since this expression behaves as 



y(x) ~ cx T for large values of x 3> 1, it reproduces the leading term of the expansion (12 1, provided t = 7. Therefore 
this form also suggests itself for a global fit of the results. 

We first restrict ourselves to the range -1< x < 1. The fit to the data for H = 1.0 X 1CT 4 MeV 5/2 , H = 1.0 X 10~ 3 
MeV 5 / 2 , and H = 1.0 X 10~ 2 MeV 5 / 2 is shown in Fig. [I] The field H varies over two orders of magnitude, yet the 
rcscalcd results fall neatly onto a single scaling curve. For these extremely small values of H, corrections to the leading 
scaling behavior are very small. In the plot, for x < 1 the differences between the scaled results are smaller than the 
width of the symbols. 

For a global fit, we use all points in the range — 1 < x < 5.3 x 10 4 . The result is dominated by the behavior at 
asymptotically large values of x. Using the scaling laws, we find from our values for the critical exponents (3, S and v 
for the critical exponent 7 

7= (2- n)v= 1.611(12), 7 = 13(5- 1) = 1.601(24), (47) 
which is still compatible with the result expected from the asymptotic behavior when we identify 

7 = t= 1.594(1). (48) 



Below we also perform a fit to the exact expression of Griffiths' expansion, including the next correction term. This 
leads to the value 7 = 1.5997(1), in perfect agreement with the result from the scaling laws. 
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FIG. 3: Comparison of the results for the Widom- Griffiths scaling function y(x), where x = t/M 1 ^ 13 and y = h/M s (M > 0). 
Shown is the result from our fit to the RG results for H = 1.0 X 10" 4 MeV 5/2 , 1.0 x 10" 3 MeV 5/2 , and 1.0 x 10" 2 MeV 5/2 
for a;-values x < 1 (black solid line, small-a; fit) and for all x values (up to approximately x — 5.3 x 10 4 ) (blue dashed line, 
asymptotic fit). There is little difference between these fits for x < 0, but for asymptotically large x- values the difference 
between these curves due to the difference in the exponents becomes quite large. The fit to the 0(4) lattice data from [57j (red 
dashed line, lattice data fit) differs from our results for large x-values due to the different values for the critical exponents. The 
result from the e-expansion [25j (green dot-dash line) is valid only for x < and displays a significantly smaller curvature than 
both our RG result and the 0(4) lattice simulation result. 



The values for the fit parameters for both the range —1 < x < 1 and the global range — 1 < x < 5.3 x 10 4 are 
given in Tab. |IV| The quality of the fit for small x is better, and it appears that the correction terms to the leading 
behavior in Griffiths' expansion for large x become important. Incidentally, the fits to the scaling function for large 
and small ir-valucs differ only in the value of the exponent r, and have the same value for the coefficient c. 

In Fig. [3] our results for the scaling function for both small and large x are compared to the fit to lattice Monte- 
Carlo results from [27] and to the result from e-expansion [25] over the interval — 1 < x < 1. For x < 0, our fits 
and the results from the lattice simulations are very close together and have a significantly larger curvature than 
the e-expansion result. For x > 0, all curves diverge quickly: The asymptotic behavior is determined by the value 
of the critical exponent 7, and the values in our calculation and the lattice simulation differ. The result |25] from 
the e-expansion satisfies the asymptotic behavior only order by order in e, but not explicitely, and is not a good 
description of the scaling function for x > 0. 



We adapt the interpolation idea of eq. (24 1 from [27 to parameterize the scaling function y(x) over the full range 
of a;- values by combining the two-parameter fit for x < 1 and the two-parameter fit for large a:-values which captures 
the asymptotic behavior determined by the exponent 7. We find that the result is not very sensitive to the exact 
point xq at which we switch from one functional from to the other, as long as < xq < 10, but smaller values give a 
slightly better interpolation. A suitable parameterization is 



"W = (i + i^a + ^ c - (1 + l) ' , + (i + i^a + ^ c ' (1+l) " with 

c s = 0.9928, t s = 1.6712, q = 0.9928, r x = 1.5941, x = 0. (49) 



We will use this expression together with the relations (111 to obtain the scaling function in the form f(z) for later 
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lo h h X 2 # d.o.f. x 2 /#d-o.f. 

0.980(5) 0.506(16) - 2.7 x 10 6 282 9586 

1.0058(6) 0.642(3) 0.135(2) 39963 281 142.3 

1.11142 0.273902 -0.015625 D. J. Wallace and R. K. P. Zia |26] 



TABLE V: Fit of the form y(x) = (1 + x) [lo + h log(l + x) + h log 2 (l + x)] . The ansatz proves to be not suitable to fit our 
results over the range — 1 < x < 1. 



as 6i X 2 # d.o.f. xV(# d.o.f.) 

0.9661(13) - - 2.3 x 10 6 283 8153 

1.2619(52) -0.2601(45) - 5.0 x 10 4 282 178.1 

1.3457(8) -0.3977(12) 0.0561(5) 277.6 281 0.9879 



TABLE VI: Coefficients for the fit y{x) = a 3 (l + x) 3/2 + 6i(l + x) + Oi(l + x) 1/2 to the 284 points in the region -1< x < 1. 
The fit with 3 parameters provides a reasonable description of our results for H — 1.0 x 10 -4 MeV 5 . 



comparison with additional results. 



2. Scaling behavior for small x 

For a determination of the scaling function for small x, we will use only the results for if = 1.0 X 1CP 4 MeV 5 / 2 in 
the range — 1 < x < 1. The corrections to the scaling behavior are small only if the external symmetry-breaking field 
is small, and this value is the smallest one for which we calculated results over a reasonably large range of values for 
the scaling variable z. 

Here we estimate the scaling corrections by comparing to the results with H — 2.0 x 10~ 4 MeV 5 / 2 . Due to the 
normalization y(0) = 1, the error is minimal at x = 0. It is largest close to x = — 1. The relative error due to scaling 
corrections is less than 0.04 in the limit x — > — 1, and for x > —0.75 it drops to less than 0.01. In the region < x < 1, 
the relative error is less than 0.0005. 

In the absence of analytical results, different ansatze for the scaling function are possible. The results [551 I2HJ US!] 
from the e-expansion suggest from the expansion around x = — 1 to attempt a fit of the form 

y(x) = (1+x) [l + h log(l + x) + l 2 log 2 (l +x) + ...], (50) 



which is only expected to hold for (1 + x) <^ 1, or from the proposed expansion (19 1 

y(x) = a 3 (l + xf' 2 + 6i(l + x) + ai(l + x) 1/2 + . . . . (51) 
Assuming that the e-expansion can be resummed, a third possibility is our phenomenological expression 

y(x) = c(l + x) T {l + di(l + x) 1/2 + d 2 (l + x) + ...). (52) 

We will take these possible expansions in turn. The results for the logarithmic expression (50 1 are given in Tab.|V] 
We find that this expression is not suitable to fit our results over the range we consider here. This is not unexpected, 



since the expansion is only valid close to x = —1. Both the expression (51 1 a nd (|52[ ) describe the results about equally 
well, the first one perhaps a little better. The coefficients are given in Tables |VI| and VII For practical purposes, both 
descriptions are equivalent. 



c r 


di d 2 x 2 


# d.o.f. 


X7(# d.o.f.) 


1.0031(5) 1.6472(20) 


3.0 x 10 4 


282 


106.3 


1.2242(27) 1.7542(14) 


-0.179(2) - 431.4 


281 


1.536 


1.390(9) 1.797(2) 


-0.332(7) 0.055(3) 77.18 


280 


0.2756 



TABLE VII: Coefficients for the fit y(x) = c(l + x) T (l + d x (l + x) 1/2 + d 2 (l + x)) to the 284 points in the region -1 < x < 1. 
Both the 3- and the 4-parameter fits are a reasonable description of the results for H — 1.0 x 10 -4 MeV 5 / 2 . 
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C2 Cl + ^3 (h 

0.681(4) 0.391(8) -0.074(4) our result 

0.674(8) 0.345(12) -0.023(5) J. Engels and T. Mendes [HI [27] 
0.746(3) 0.19(1) 0.061(8) A. Cucchieri and T. Mendes [28] 
0.530 0.528 - D. J. Wallace and R. K. P. Zia [H] 

TABLE VIII: Fit to the equation of state in inverted form, 1+x — c 2 y 1 '' 2 + (ci+d 3 )y+d2y 3 ^ 2 , from our results for H — l.Ox 10 -4 
MeV 5/2 in the region < y < 1.5, compared results from the e-expansion and from 0(4) spin model lattice calculations. 



7 ci £2 x 2 # d.o.f. xV# d.o.f. 

1.5948(4) 0.9855(25) - 2144 204 10.51 

1.59973(4) 0.94294(28) 1.196(8) 4.186 203 0.02062 



TABLE IX: Fit to Griffiths' expansion y(x) = x 1 (a + c 2 x~ 213 + c 3 x~ 4/3 + ...), for x > 0, M > M > 0. We fit only to the 
results for H = 1.0 X 10" 4 MeV 5/2 and only to points with x > 100. 



In [27] and |26j . the equation of state is parameterized in inverted form as x = x(y) 7 and the terms expected from 
the e-expansion in d = 3 (e = 1) are 

l + x = c 2 y 112 + (ci + d 3 )y + d 2 y 3/2 + . . . . (53) 

To make a direct comparison to our results easier, we also used this form of the equation for a fit to the results for 
if = 1.0 x 1(T 4 MeV 5 / 2 in the range 0<v < 1.5. 



The coefficients are given in Table VIII The coefficient d 2 has not been calculated in the e-expansion. For the 
comparison to the lattice results, please keep in mind that our results have an additional systematic error not included 
in the fit uncertainty given in the table. While strictly speaking the results from the spin model lattice calculation 
[181 I27j and from our results do not agree within the errors, they are remarkably consistent with each other, in 
particular when compared with the perturbative results from the e-expansion: The values for the leading coefficient, 
c 2 , agree within 1%. The values for the sub-leading coefficent, c.\ + d 3 , differ only by about 10%. In contrast, both 
our results and the lattice results differ from the perturbatively calculated coefficients from the e-expansion by more 
than 20 - 30%. 



3. Asymptotic scaling behavior 



Griffiths' expansion (12) can be used to describe the results for large values of x, and to determine the critical 
exponent 7 to high accuracy. A comparison to the value obtained with the scaling laws from the other exponents 
serves as a consistency check for our results. We use for the fit 

y(x) = cix 7 + c 2 x^' 2)3 + c^- Ap + . . . = x-' (ci + C2X~ 20 + c^x^ + . . .) (54) 

with 7 as a free parameter. In order to retain the hierarchy of the corrections in the expansion, we use the previously 
determined value j3 = 0.4030(30) and keep it fixed. Since 7 — 4/3 ^7 — 2v«0 according to the d = 3 scaling laws, 
only the exponents of the first two terms in Griffiths' expansion are positive and contribute significantly for large x. 
The additional terms with n > 3 are small corrections for large x, but diverge for x — > where the expansion is no 
longer valid. 

We fit again only to the results for H = 1.0 x 10~ 4 MeV 5 / 2 and estimate scaling corrections from a comparison 
to the results for H = 2.0 x 10 -4 MeV 5 / 2 . Only points with x > 100 are included, so we can be certain to be in a 
region where the expansion is valid. The results for the coefficients and the exponent 7 are given in Tab. |LX| The 
inclusion of more than the first two terms does not lead to a meaningful improvement of the fit quality. As already 
noted above, the results for the critical exponent 7 = 1.5997(1) is in perfect agreement with the determination from 
the scaling laws. 



B. Scaling behavior for large fields 



After having determined the scaling function from results with small values of the symmetry-breaking field H, we 
now compare these results to the scaling behavior at much larger values of H. For this purpose, we also switch to the 
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FIG. 4: The parameterization of the scaling function y(x) in terms of the scaling variables x and y pro vides an implicit 
parameterization of the scaling function f(z). Shown is a comparison of f(z) from the interpolated fit eq.(49 1 with the rescaled 
results for H — 1.0 x 10~ 4 MeV 5 / 2 (blue circles, for clarity not all points are shown). The fit is perfect on the scale of the 
plot. The asymptotic behavior of y(x) for x — > oo determines the behavior of the fit for z > 0.5. Using the parameterization 
in terms of only the small- a; behavior underestimates f(z) in this region. In this region, the results from the e-expansion show 
the largest deviation. 



more intuitive description by means of the scaling function f{z). 

With the help of Eqns. (11) and the critical exponents, the parametrization Eq. ( |49| implicitly also provides a 
parameterization of the s calin g function in the form f(z). In Fig. [4j the interpolation function f(z) obtained from 
the parameterization Eq. (49) is compared to the rescaled results for the order parameter M/h 1 / 5 for H — 1.0 x 1CP 4 
MeV 5 / 2 as a function of z. On the scale of the plot, the agreement is perfect. The scaling form f(z) is less sensitive 
to small changes in M than the Widom-Griffiths scaling form y(x) and thus the agreement appears to be even better 
than for the fit in the form y(x). 

Using the appropriate values for the critical exponents, we can also translate the scaling functions from the e- 
expansion and from the lattice simulations [37] into this form. The comparison to our results for H = 1.0 x 10~ 4 
MeV 5 / 2 is shown in Fig. [H]for the e-expansion and in Fig. [6] for the lattice results. Since f(z) — ► (—z)P for z — ► — oo, 
the differences in the values of the critical exponents lead to different behavior of the scaling functions. To illustrate 
this, the asymptotic behavior for the different values of /3 is also shown in the plots, and the scaling functions can be 
seen to approach the asymptotic functions. 

The result from the e-expansion shows the largest deviation from our results for z > 0.5, which corresponds to the 
deviation observed in y(x) at large x-values, where the expansion is not a good description of the scaling function. 

The agreement between our results and those from the lattice simulations is much better. Most of the difference 
at large negative z-values can be attributed to the different values of the critical exponents. Overall we consider the 
agreement quite satisfactory. 

The scaling behavior of the order parameter can be demonstrated very nicely in this scaling form with our results. 
Having first confirmed the agreement between our scaling function f(z) and the rescaled results for H = 1.0 x 10~ 4 
MeV 5 / 2 in Fig. [l] we now compare these results to results for a wide range of values for H. 

In Fig. [7] we consider small values of H in the range from H = 1.0 x 10~ 4 to H = 1.0 x 10~ 3 MeV 5 / 2 . In the 
left-hand panel, the order parameter M is plotted as a function of the reduced temperature t for the different values 
of H. In the right-hand panel, the rescaled order parameter M/h 1 / 5 is plotted as a function of z = t/h 1 '^ for the 
same H values. After rcscaling, the curves for all values of H fall onto a single line and are indistinguishable at the 
scale of the plot. This agrees with our observation in the Wisdom-Griffiths scaling analysis, where scaling corrections 
were also negligible in this 7J-range. 

However, for larger values of H, corrections to the scaling behavior soon become apparent. In Fig. [8] we show 
results in the range from H = 1.0 MeV 5 / 2 to H = 1.0 x 10 3 MeV 5 / 2 . A gain both the order parameter as a function of 
temperature and the rescaled order parameter as a function of z are shown. For comparsion, the rescaled results for 
H = 1.0 x 10~ 4 MeV 5 / 2 and the asymptote (— z) 13 are shown with the rescaled results. While the curves for different 
values of H still collapse after rescaling in the vicinity of the critical temperature, the deviations from the scaling 
function become quite large, and they are already plainly visible for fields of the order H = 10.0 MeV 5 / 2 . 

This result is actually not very surprising: For small values of the field H, only a small region in t around the critical 
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FIG. 5: Plot of the scaling function f(z) as a function of z. Shown are results from our calculation for H = 1.0 X 10" 4 MeV 5/2 
(black solid line) and from the e-expansion calculation to 0(e 2 ) of Brezin et al. [55] (green dashed line). For comparison, we 
also plot the expected asymptotic behavior f(z) — > (— z) 13 for j3 = 0.390 (e-expansion) (dot-dashed green line) and for j3 = 0.403 
(our result) (dotted black line). The behavior of the scaling functions for large negative values of z is again determined by the 
different values of the critical exponent (3. The e-expansion result is expected to provide a good description only for t < 0. We 
find that indeed the deviation from our scaling function is large for z > 0.5. 
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FIG. 6: Plot of the scaling function f(z) as a function of z. Shown are results from our calculation for H — 1.0 x 10 4 
MeV 5//2 (black solid line) and from the calculation of Mendes and Engels 127, (using the parameterization as given in |18| ) (red 
dashed line). For comparison, we also plot the expected asymptotic behavior f(z) — > (— z) 13 for /3 = 0.380 (Mendes and Engels) 
(dot-dashed red line) and for (3 = 0.403 (our result) (dotted black line). As expected, the behavior of the scaling functions for 
large negative values of z is determined by the different values of the critical exponent j3. 
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FIG. 7: Order parameter as a function of the reduced temperature t for very small values of the external field H from 1.0 x 10~ 
to 1.0 x 10~ 3 MeV 5/2 (left panel), and for the rescaled order parameter M/h 1/s as a function of z = t/h 1/{J3S) (right panel). 
The i-ranges for the different values of H are chosen such that z covers the range —10 ... 10 after rescaling. On the scale of the 
plot in the left-hand panel, a deviation from the scaling behavior is not visible in this range of values. 



temperature t = is rescaled and becomes a sizable region in z. Therefore we can expect that all results are in the 
scaling region around T c and scaling deviations are small over a large region in z. For example, for H — 1.0 x 10~ 4 
MeV 5 / 2 , the rescaled results approach the expected asymptotic behavior at least down to z = —30. In contrast, for 
large values of H, a much wider region in t is mapped onto the same range of z-values, and corrections to the scaling 
behavior are bound to appear, since points far away from the critical point are included. In Fig. [5J one can see how 
the t-region has to be increased with H to keep the z-range constant after rescaling. Results from the boundary of 
the scaling region at large —t then lead to deviations from the expected asymptotic scaling behavior at large —z. 

The deviations from the leading-order scaling behavior are even more starkly visible from a plot of the results in the 
Widom-Griffiths scaling form as y(x). This is shown for values of H over 7 orders of magnitude, from H = 1.0 x 10~ 4 
MeV 5/2 to H = 1.0 x 10 3 MeV 5/2 , in Fig. [9] for the range -1< x < 1, and in Fig. [ro] for the range -1<i<5.0x 10 3 . 



On the scale of these plots, corrections become significant for fields H > 1.0 x 10^^ MeV 5 / 2 . 

In conclusion, we find that we can extract a parameterization of the scaling function in Widom-Griffiths form 
y = y(x) from our results for the order parameter at small values of the field H . The scaling function f(z) can be 
obtained from this parameterization, and the agreement with the order parameter in this scaling form is also very 
good. These results are also in satisfactory agreement with the results for the scaling function from the 0(4) spin 
model lattice simulations [27], and the principal differences for large negative values of the scaling variable z can be 
explained by the difference in the values of the critical exponents. 

Beyond small values of the external field H , we find consistently in both scaling forms, y(x) and f(z) , that corrections 
to scaling become large and the results deviate from the scaling function away from the critical temperature. The 
scaling deviations can be discerned more readily in the Widom-Griffiths scaling form. 



VI. SUSCEPTIBILITY 



The susceptibility affords us an additional opportunity to test our results for the scaling function. As outlined in 
the discussion of the scaling behavior in Sec. [TTJ the scaling function for the susceptibility is completely determined 
by the scaling function for the order parameter and by the critical exponents: 



1 



or [H M 5 - l x] =Sy(x)--xy'(x) 



From a view focussing on the effective potential, this relation is a priori far from obvious: the order parameter is 
determined as the minimum from the first derivative of the effective potential, the susceptibility from the second 
derivative, which is related to the masses of the longitudinal fluctuations and the four-point interaction. For this 
reason we consider this to be a significant test not only of our result for the scaling function, but also of our approach 
in general. 
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FIG. 8: Results for the order parameter as a function of the reduced temperature t for different values of the field H (left-hand 
column), and for the rescaled order parameter M/h 1 / 1 as a function of the scaling variable z = t/h 1 '^^ (right-hand column). 
For comparison, the scaling function obtained from the result for H = 1.0 x 10~ 4 MeV 5//2 (solid black line) and the asymptotic 
behavior of the scaling function {—z)^ (dashed black line) are shown with the rescaled results. Note that the axes on the figures 
for the rescaled order parameter in the right-hand column all have the same scales, while the scales of the axes on the plots 
with the unsealed results differ widely for the different -ranges. 
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FIG. 9: Comparison of scaling results over a wide range of values (7 orders of magnitude) for the symmetry-breaking field H . 
On this scale, deviations from the scaling behavior become readily apparent for fields of the order of H = 1.0 x 10 _1 MeV 5 . 
For comparison, the fit y(x) = c(l + x) T with c = 0.9228(24) and r = 1.6712(27) to the results for H = 1.0 x 10 -4 , 1.0 x 10" 3 , 
and 1.0 x 10 -2 MeV 5//2 for x < 1 is also shown (solid black line). For clarity, not all points are shown. 



a;-range 


c t x 2 


# d.o.f. 


*7# d.oi. 


-1 < x < 1 


1.0009(19) 1.7018(33) 371.1 


838 


0.4428 


-1< x < 5.3 x 10 4 


0.9802(27) 1.5954(12) 3754 


1551 


2.420 



TABLE X: Fit to the inverse rescaled susceptibility for H = 1.0 X 10~ 4 , 1.0 x 10" 3 , 1.0 x 10" 2 MeV 5/2 with ansatz for 
the equation of state of the form y(x) — c(l + x) T . The scaling corrections are estimated from the spread of the rescaled 
susceptibility. Shown are separate fits for the range — 1 < x < 1 and for — 1 < x < 5.3 x 10 4 . 

In the following, we will proceed in the same way as for the order parameter. We will first determine the scaling 
function for the susceptibility in Widom-Griffiths form from small values of the eternal field H, and then proceed to 
larger values. We will compare the determination of the scaling function from the susceptibility to the determination 
from the order parameter. For completeness, we also compare to the scaling function obtained from the lattice scaling 
function for the order parameter. 

A. Scaling function in Widom-Griffiths parameterization 

1. Phenomenological ansatz for the scaling function 



As for the order parameter, we fit a two-parameter ansatz y(x) — c(l + x) T for the equation of state to the results 
for if = 1.0 X 10~ 4 MeV 5/2 , if = 1.0 X 10~ 3 MeV 5/2 , and H = 1.0 X 10~ 2 MeV 5/2 . The size of the scaling corrections 
is estimated from the spread of the rescaled results for these values of H . Inserting the ansatz in the expression for 
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FIG. 10: Comparison of scaling results over a wide range of values (7 orders of magnitude) for the symmetry-breaking field 

H. Deviations from the scaling behavior are clearly visible for fields of the order of H = 1.0 x 10 _1 MeV 5 / 2 . For comparison, 
the asymptotic fit y(x) = c(l + x) T with c = 0.9228(24) and r = 1.5941(8) to the results for H = 1.0 X 10~ 4 , 1.0 x 10" 3 , and 

I. x 10~ 2 MeV 5//2 is also shown (dashed blue line). For clarity, not all points are shown. 



the rescaled susceptibility, one obtains 

Sy(x) - ~xy'(x) = c(l + xY (s - ~ , (55) 

which we use to fit the inverse rescaled susceptibility. The critical exponents f3 and 5 are kept fixed at the previously 
determined values. 

The results for the separate fits in the range — 1 < x < 1 and in the range — 1 < x < 5.3 x 10 4 are shown in Fig 
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and Fig. [12j Scaling for the susceptibility, here over two orders of magnitude for the field H , is clearly observed, with 
small corrections at these values of H . By comparing the results at large x in Fig. [12] for the susceptibility to the 
results in Fig. |2j for the order parameter, one can already see an indication that Sy(x) — jjxy(x)' ~ y(x) as expected 
from Griffiths' expansion. We will check this more rigorously below when we analyze the asymptotic behavior. 
The coefficients for the two-parameter fits are shown in Table [X] They can be compared directly to the coefficients 



from the fit to the rescaled order parameter in Table IV While the agreement is perfect (better than 0.1%) only for 
the exponent r for the large-x region, which corresponds to the critical exponent 7, all parameters agree to better 
than 2%. 

We can find a global expression for the rescaled susceptibility that interpolates between the small- a; and the large- a; 
behavior, just as for the order parameter. We first form the combinations 

Sy(x) - ^xy'(x) 
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FIG. 11: Results for the inverse scaled susceptibility [H M s - 1 x] 1 for H = 1.0 x 1(T 4 , 1.0 x 10" 3 , and 1.0 x 10" 2 MeV 5/2 
(not all points are shown). We estimate an error due to the scaling violations from the spread of the results. Shown are in 
addition the fits for small x- values x < 1 (solid black line), and from all values up to approximately x = 5.3 x 10 4 (dashed blue 
line, asymptotic fit). 
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TABLE XI: Coefficients for the phenomenological ansatz for the interpolated equation of state. Values in the upper row are 
obtained from a fit to the susceptibility \, values in the lower row from the order parameter M. 



for the large- a; and the small- x region separately, and then combine them with interpolation factors: 

\M & - l H xY l = r, { \t X °n to c s (l + *r-j(l v . 

L J (1 + x a ) 2 + (1 + x) 2 v ' V /3Sl + x / 



, cAl + xY'S 1- ~ . (56) 

(l + x ) 2 + (l + a;) 2 v ' V (361 + x J 

For the coefficients, we give both the values obtained directly from the susceptibility, and the values obtained from 
the order parameter in Tab. |XI| For all practical purposes, these functions are indistinguishable from one another. 
This is illustrated in Fig. |13| Thus the predictions of scaling are borne out by our results, and the scaling function 
for the susceptibility is already predicted by the scaling function for the order parameter. 

This result can be compared to the interpolated parameterization of the lattice results from Engels and Mendes. 
We use the parameterization x(y) for small and large x-values obtained from the order parameter. For both regions 
we first form the combinations 



P x'{y) ' 



where in this case y is treated as the independent variable. We then use the interpolation prescription eq. (24 1 to 
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FIG. 12: Results for the inverse scaled susceptibility [H M s - 1 x] 1 for H = 1.0 x 10~ 4 , 1.0 x 10" 3 , and 1.0 x 10" 2 MeV 5/2 
(not all points are shown). The fit to all data points is also shown (asymptotic fit, blue dashed line). From comparison to 
the equation of state y(x) in Fig. |2jit is apparent that for large x indeed 6y(x) — ^xy(x)' ~ y(x) as expected from Griffiths' 
expansion and the scaling laws (see discussion of the asymptotic behavior below). 



c r 


d\ d 2 x 2 


#d 


o.f. X 2 /# d.o.f. 


1.00324(26) 1.6546(17) 


3.0 x 10 4 


282 


109.2 


1.299(10) 1.7886(43) 


-0.229(4) - 2762 


281 


9.827 


1.5344(49) 1.8472(13) 


-0.428(4) 0.081(2) 98.51 


280 


0.3518 



TABLE XII: Coefficients for the fit to the inverse scaled susceptibility in the region —1 < x < 1, using the ansatz y(x) = 
c(l + x) T (l + di(l + x) 1/2 + d 2 (l + x)) for the equation of state. 



obtain a result valid over a wide y-range. For the comparison in terms of x and y and for small x-values, the resulting 
curve is also shown in Fig. [13] as the lattice data fit. The difference in the slope is to a large part once again due to 
the difference in the values for the critical exponents. 



2. Scaling behavior for small x 

We will again determine the scaling function for small x as accurate as possible from the results for H = 1.0 x 10 -4 
MeV 5 / 2 in the range — 1 < x < 1, where scaling corrections are small. We estimate the scaling corrections for the fit 
to the smallest-field results H = 1.0 x 10~ 4 MeV 5 / 2 by comparing the rescaled results to those of with H = 2.0 x 10~ 4 
MeW 2 . Between —0.5 < x < 1, the scaling violation is less than 0.01%, and becomes sizable (~ 5%) only close to 
x = — 1. Using these errors, we can judge how well the leading-order scaling behavior fits the results. 

We use the same ansatze as for the order parameter. With the ansatz y(x) = c(l + x) T (l + di(l + x) 1 / 2 + ^(l + x)) 
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FIG. 13: Comparison of different fits for the universal function 5y(x) — hxy'(x) which describes the inverse rescaled sus- 
ceptibility. We show the function with parameters obtained directly from a fit to our results for the rescaled susceptibility 
(blue dot-dashed line), from our results for the order parameter (black solid line), and from the order parameter in the lattice 
simulation [27] (note that the values for the critical exponents are different). For all practical purposes, there is complete 
agreement between our curves with parameters from the order parameter and with those from the susceptibility. 



«3 02 gi X 2 # d.o.f. x7(# d.o.f.) 

0.98125(49) - - 3.4 x 10 6 283 12100 

1.2855(17) -0.2839(16) - 7746 282 27.47 

1.3678(16) -0.4409(29) 0.0751(15) 196.9 281 0.7007 



TABLE XIII: Coefficients for the fit to the inverse scaled susceptibility in the region —1 < x < 1, using the ansatz y(x) = 
03(1 + t) 3 / 2 + 02(1 + x) + ai(l + x) 1 ^ 2 for the equation of state. 



for the equation of state, one finds for the scaling function of the susceptibility 



1 



Sy(x)- -xy'(x) = S c(l + x) T 1 + d^l + x) 1/2 + d 2 (l + x) 



-5c(l + x) T 



138 (1 + x) 



t + (t+ -)di(l + x) 1 ' 2 + (r + l)da(l + x) 



(57) 



We take the critical exponents from our original determination, and thus fit the result in this functional from with the 
same number of parameters as for the order parameter. The results are given in Table |Xn] They can be compared 
directly to the results from the order parameter in Tab. |VII| The agreement with the parameters determined from 
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7 ci C2 c 3 X 2 # d.o.f. xV# d.o.f. % - l/(ffi)) 

1.5844(12) 1.081(9) - - 9918 204 48.62 1.0411 

1.59988(5) 0.94158(42) 1.2794(34) - 4.2109 203 0.02074 1.00277 



TABLE XIV: Coefficients for the fit to the inverse susceptibility in the region x > 100, using the asymptotic expression 
y(x) — x 1 '(ci + C2X~ 213 + czx~ il3 + . . .) for the equation of state. The deviationfrom 1 in the last column indicates the violation 
of the scaling law for the critical exponents with the value of 7 from the fit. 



the order parameter is best for the fit with only two parameters: we find c = 1.0032(3) and r = 1.655(2) from the 
susceptibility, and c = 1.0031(5) and r = 1.682(3) from the order parameter. The agreement for the coefficient c is 
perfect, and the values for the exponent r agree within 2%. Even though the agreement of the coefficients becomes 
worse when additional corrections are included, the functions with parameters from the order parameter and from 
the susceptibility remain within 0.4% of each other over the interval —1 < x < 1. 

The ansatz y(x) = a^(\ + x) z / 2 + 02(1 + x) + ai(l + x) 1 / 2 for the equation of state leads to the corresponding 
expression 
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1 1 



(58) 



The coefficients are given in Table XIII They should be compared to those from the fit to the order parameter 
in Table VI Again, the leading-order coefficients agree within 2% in the one-parameter, and within 5% in the two- 
parameter fit. While the agreement between the coefficients is not perfect, the agreement between the scaling functions 
as such is actually better: For the three-parameter fit, the scaling functions fitted directly to the susceptibility and 
the one calculated from the fit to the order parameter are within 0.01 of one another over most of the interval 
(—0.9 < x < 1), i.e. the relative errors arc less than 0.4% for most of the interval (—0.9 < x < 1). Since the function 
itself approaches for x — > — 1, the relative errors become large close to x = — 1. 

Despite the differences in the coefficients that appear when additional corrections are included, the determination 
of the scaling function from the susceptibility and the one from the order parameter are completely equivalent for all 
practical purposes. 



3. Scaling behavior for asymptotically large x 
Using Griffiths' expansion eq. (|12[) in the scaling function of the susceptibility, one finds that the leading order 



behavior of the scaling functions for the order parameter and the susceptibility agree for large x (see eq. (14)), 
provided the scaling law 7 = 0(5 — 1) holds. In the following, we test this by fitting Griffiths' expansion to the 
results from H=l.0x 10~ 4 for x > 100, once again with the deviation from the scaling behavior estimated from the 
difference to the results for H = 2.0 x 10~ 4 and used to assess the quality of the fit. As for the order parameter, we 
keep the values of the critical exponents and 5 fixed to the previously determined values, and treat 7 as a free fit 
parameter. 

The coefficients of the fit are given in Table |XIV| They should be compared to the results for the coefficients of the 
equation of state in Tab. |IX| Because of the scaling relations between the critical exponents, the leading-order terms 
of the scaling form for the order parameter, y(x) ~ cicc 7 , and for the inverse susceptibility, 5y(x) — jjxy'(x) ~ cix 1 
are expected to coincide exactly. 

With two fit parameters, the fit results for the exponent 7 agree to within 0.0004, and the coefficients c\ to within 
9%. Part of the difference for c\ is explained by the fact that the exponent 7 satisfies the scaling law only to within 
4%, which translates directly into a deviation for the coefficient c\ and explains about half of the difference. 

For the three-parameter fit (with 7, c\ and ci as fit parameters), the agreement between the equation of state y(x) 
from the order parameter and from the susceptibility is perfect, the values for the critical exponent 7 agree to within 
0.0002, and the scaling law 7 = 0(5 — 1) is satisfied to better than 0.3% (see deviation from 1 in the last column of 
the table). 

We find the degree to which the relations between the scaling functions are satisfied in our calculation truly 
remarkable. On the one hand, we calculate the order parameter from the minimum of the effective potential. On the 
other hand, the longitudinal susceptibility x is calculated from the mass m a of the longitudinal fluctuations, which 
involves the curvature of the effective potential and is a priori independent of the order parameter. Nevertheless, 
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FIG. 14: The parameterization y — y(x) of the equation of state in terms of the scaling variables x and y also provides 
implicitly a parameterization of the scaled susceptibility Hoh 1-1 ^ 6 % = § \j{ z ) — jj m terms of the scaling variable z. The 
scaling function (black solid line) with parameters obtained from the order parameter, not the susceptibility itself, is compared 
to the results for H = 1.0 x 10~ 4 MeV 5//2 (blue circles, for clarity not all points are shown). Apart from the slight difference 
at the peak the fit is almost perfect. Within the linewidth, the fit with parameters obtained from the order parameter and the 
one with parameters obtained from the susceptibility directly (not shown) are indistinguishable in this plot. 



we find that the values for the critical exponents and for the leading coefficients in the asymptotic expansion of the 
equation of state agree perfectly. 

In this region, this relation has been tested neither by results from the e-expansion, nor by results from lattice 
simulations. For the e-expansion, the large-x region is not accessible, since the asymptotic behavior cannot be 
calculated explicitly, but only order by order in e. In spin model lattice calculations, the equation of state has also 
not been calculated this far in x and the asymptotic behavior has not been tested to this extent. Engels and Mendes 
[27) confirm the large-?/ behavior in the scaling form x(y) up to x < 150, after inverting the equation of state y(x) 
approximately. 



B. Scaling for large fields 



We have determined the scaling function for the susceptibility in the previous sections from results at small values of 
H, and we were able to confirm the expected relations with the scaling function for the order parameter. Proceeding 
as for the order parameter, we now turn to larger values of H where scaling corrections become important. 

The representation of the results in terms of Griffiths' scaling variables x and y is not very intuitive, and we 
therefore translate the results back into the scaling variable z and the scaling function f(z). 



In Fig. 14 



susceptibility H h x -V*x = \ /(*) - ff'(z) 



for H = 1.0 x 10~ 4 MeV 5 / 2 is shown as a function of 



the rescaled 
together with 

the fit with parameters obtained from the order parameter. Apart from the slight difference at the peak, which is 
most sensitive to scaling corrections, the scaling function agrees very well with the results for the susceptibility. This 
really confirms a prediction, since all parameters for the curve are already determined from the order parameter, and 
there is no freedom left. A direct fit to the susceptibility in a second step only confirms how well the parameters agree 
with each other. 

A comparison of our results for H = 1.0 x 10~ 4 MeV 5 / 2 to the scaling function for the susceptibility predicted from 
the lattice spin model is shown in Fig. 
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For z < (in the phase with large symmetry breaking), the agreement is 
very good, but in the phase with largely restored symmetry, the susceptibility predicted by the order parameter lattice 
scaling function is larger than the one we calculate in the RG. Overall, we consider the agreement still satisfactory. 
We now demonstrate the scaling behavior of the suscep tibi lity by plotting the susceptibility rescaled in the form 



H h' 
H 



" 1/,<5 X as a function of the scaling variable z. In Fig. 16 results for small values of H = 1.0 x 10~ 4 MeV 5 / 2 to 
1.0 x 10~ 3 MeV 5 / 2 are shown. In the left panel, the susceptibility as a function of the reduced temperature t is 
plotted for different values of H. In the right panel, the rescaled susceptibility Hoh 1 ~ 1 / 5 x IS plotted as a function of 
z = i/ftVW for the same values of H. As for the order parameter, we observe that the results collapse onto a single 
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z = t/ti 



FIG. 15: Comparison of the parameterization of the scaling function from the lattice results !27] (red dashed line) to our 
results for H = 1.0 x 10~ 4 MeV 5 ^ 2 (black solid line). The agreement is very good below the critical temperature, but deviations 
appear above T c . 




FIG. 16: Susceptibility as a function of the reduced temperature t for very small values of the external field H = 1.0 x 10 
MeV 5/2 to H = 1.0 x 10" 3 MeV 5/2 (left panel), and rescaled susceptibility iW 1/<5 X as a function of z = t/h 1/(0S) (right 
panel). The f-ranges for the different values of H are chosen such that z covers the range —10 ... 10 after rescaling. Within the 
width of the symbols, the rescaled curves coincide exactly. Deviations from scaling are not visible at this scale. 
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FIG. 17: Results for the longitudinal susceptibility \ as a function of the reduced temperature t (left-hand column), and for 
the rescaled susceptibility Hoh 1 ~ 1 ^ s x as a function of z = t/h 1 '^ 6 ' (right-hand column) for different values of the field H . 
For comparison, the scaling function obtained from the results for H = 1.0 x 10~ 4 MeV 5 ^ 2 is plotted with the rescaled results 
(black dashed line). Note that the axes for the rescaled results in the left column are kept at the same scale, while the ranges 
in t for the unsealed results are very different for different i!f -ranges. 
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i — i — i — i — | — i — i — i — r 




FIG. 18: Comparison of the scaling results for the inverse susceptibility [HqM 1 \] over a wid e range of values for the 
symmetry-breaking field H (7 orders of magnitude). The result of a fit with the ansatz eq. (561 to the susceptibility for 
H = 1.0 X 1(T 4 MeV 5/2 , H = 1.0 X 10" 3 MeV 5/2 , and H = 1.0 X 10" 2 MeV 5/2 is shown as a solid line, 
of H, the deviations from the leading-order scaling behavior are apparent. 



For the larger values 



curve after rescaling. No scaling corrections are discernible for these small values of H. 

However, as we increase the values for H, corrections to scaling soon become apparent. In Fig. |17| we show our 
results for the range H = 1.0 MeV 5 / 2 to H = 1.0 x 10 3 MeV 5 / 2 . Both the susceptibility as a function of t and the 
rescaled susceptibility as a function of z are shown. For comparison, we also plot the rescaled results for H = 1.0 x 1CP 4 
MeV 5 / 2 . Deviations from the scaling behavior are already significant for H = 10.0 MeV 5 / 2 , in agreement with our 
observations for the order parameter. 

The corrections to the scaling behavior become more obvious when we again plot the results in Widom-Griffiths 



scaling form. In Fig. 18 the rescaled results for H = 1.0 X 10 4 to H = 1.0 x 10 3 MeV 5 / 2 (H covers 7 orders of 
magnitude) are shown as a function of x for the range — 1 < x < 1. In this plot, the corrections already appear 
significant for H > 1.0 X 10~ 2 MeV 5/2 . 

In conclusion, we were able to obtain the scaling function for the susceptibility in Widom-Griffiths scaling form from 
the RG results at small values of H . This scaling function satisfies the expected relations with the scaling function 
obtained from the order parameter. In fact, we find that the scaling behavior of the susceptibility is described perfectly 
by the equation of state obtained from the order parameter. The scaling behavior as a function of the scaling variable 
z can be obtained from the parameterization as well. Beyond small values of the field H, we find in both scaling forms 
that corrections to scaling become large and that the deviations from the scaling functions can become significant. 
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VII. MASSES 

So far we have discussed the universal scaling behavior which is expected to apply to all systems with a critical 
point governed by 0(4) symmetry But for any scaling analysis, it is an important question how large the scaling 
region around the critical point actually is. This is not a universal property, but very much dependent on the system 
in question. Because the normalization constants To and Hq are not universal but depend on the details of the system, 
a direct extrapolation of the scaling region from one system to a different one is not possible. 

Universality only applies to the long-range behavior on small momentum scales, but not to the details of the short- 
range interactions. For different systems within the same universality class, on short scales in the UV regime the 
physics can be very different, and thus it is also not possible to fix scales for the scaling region from the UV behavior. 

However, in the context of lattice simulations, a very important technique for the scaling analysis arc finite-size 
studies. For a finite-size analysis, the length scale L provides an additional point to compare long-range behavior of 
different systems, where the systems are comparable. The quantity of interest for this behavior is the ratio of the 
correlation length to the system size, £/L. Since the mass of the longitudinal fluctuations is bounded by the mass of 
the pseudo-Goldstone bosons, m a > m^, a useful measure is given by m^L. 

We have studied finite-size scaling for the 0(4) model with the same choice of parameters as used in this work 
122] , Using the scale set by our choice of A, we find that the fields required to map out the finite-size scaling 
behavior for volume sizes of a few fm are much larger than the ones employed in this calculation. Consequently, we 
observe large corrections to scaling in these calculations. 

For the current choice of parameters, the masses are shown in Fig. |19| Compared to absolute values on a hadronic 
scale, the masses of the fluctuations in the scaling region in this calculation are very small, and scaling corrections 
quickly become large for larger masses. However, it remains possible that we can reach larger, more realistic pion 
masses while still remaining inside the scaling region. A change in the scale A changes the critical temperature and 
the values of the normalization constants T and H . This could increase the size of the scaling region, while the 
absolute value of and the dimcnsionlcss product m^L could be kept constant. 

The consequences of these observations for the analysis of QCD lattice data are less clear. Since direct comparisons 
are inadvisable, conclusions must remain somewhat speculative. We find in our results for the order parameter for 
large values of H that the results for the susceptibility as well as the order parameter still appear to scale, i.e. they 
are close together after rescaling, but show a large deviation from the scaling function. Similar behavior is observed 
in some lattice simulation studies, where scaling with the critical exponents seems to take place for the peaks in the 
chiral susceptibility, but no agreement with the 0(4) scaling functions is found [StJl I5TI 152] . 

VIII. CONCLUSIONS 

The 0(4) scaling function in three dimensions is important for the scaling analysis of systems in this universality 
class. QCD with two flavors is expected to fall into this class, if the phase transition is second order for two massless 
quark flavors. In QCD lattice simulations, the symmetry is broken by frequently large quark masses and a scaling 
analysis is necessary to determine the actual order of the phase transition. In addition, since lattice simulations are 
performed in a finite volume, finite-size scaling analysis is an important tool. Reliable knowledge about the scaling 
functions improves the power of this analysis. 

We have investigated scaling in the 0(4) model in d = 3 with a non-perturbative Renormalization Group (RG) 
calculation. In contrast to many earlier investigations of scaling with functional RG methods, we have explicitly 
included an external symmetry-breaking field H. Due to the presence of this field, we chose not to work in a scale-free 
formulation, but retain the dimensions of all quantities. The scale for the calculation is set by our choice A = 1.0 
GeV for the initial RG scale. 

We work throughout in a local potential approximation in which the anomalous dimension vanishes, r\ = 0. To 
ensure consistency of the analysis, we determine the critical exponents f3, S, and v, and find good agreement with other 
RG calculations in this approximation. The values differ systematically from those of lattice Monte-Carlo simulations 
and from Functional RG calculations with momentum-dependent couplings, which can be explained by the vanishing 
anomalous dimension. 

We determine the scaling function for the order parameter in Widom-Grifhths scaling form y(x) and as a function 
f(z) of the scaling variable z = t/h 1 ^^'. We find good agreement with the results from the 0(4) spin model lattice 
Monte-Carlo simulations of Engels and Mendes. Differences of the scaling function for the asymptotic behavior are 
explained by the difference in the values for the critical exponents. In addition, we explicitly check Griffiths' expansion 
of the equation of state for asymptotically large values of the scaling variable, and we recover the critical exponent 
7 from the asymptotic behavior. The value found in this way is in complete agreement with the one obtained with 
the scaling laws. We obtain a parameterization of the scaling function in Widom-Griffiths form y(x) which provides 
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FIG. 19: Mass of the transversal fluctuations (pseudo-Goldstone bosons) vs. reduced temperature t (left column), and 
mass of the longitudinal fluctuations m CT vs. t (right column) for different values of the symmetry-breaking field H. For our 
0(4)-model with a scale set by the UV cutoff A = 1.000 GeV, large scaling violations appear already for values of the external 
filed H which correspond to m v of less than 10 MeV. 
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a very good description of our scaling results for small values of the external symmetry-breaking field H. This 
parameterization, together with the critical exponents, also provides a parameterization of the scaling function f(z). 
This result can be used for comparison in a scaling analysis, although come caution is warranted due to the systematic 
error in the values of the critical exponents, which is also reflected in the scaling function. 

For large values of the symmetry-breaking field, we still observe scaling, but scaling corrections quickly become 
large. In terms of the scaling variable z, the region in which the scaled results fall onto the scaling curve shrinks 
considerably with increasing H . We have covered 7 orders of magnitude in H 7 from perfect scaling behavior to a region 
where scaling violations become quite large. We observe that the Widom-Griffiths scaling form is more sensitive to 
scaling corrections than a rescaling as a function of z. 

The scaling function for the longitudinal susceptibility is obtained from a direct calculation of this susceptibility. 
We confirm an important consequence expected from the critical scaling behavior and show that this scaling function 
is already given by the scaling function of the order parameter. We find very good agreement between the parameters 
for the equation of state y = y(x) obtained from fits to the order parameter and the susceptibility. We further confirm 
that the leading large- a; behavior of the scaling functions for the order parameter and the susceptibility coincide, 
which is predicted by the scaling relations and Griffiths' expansion. 

Using the scaling form /(z) for the scaling function, we find that the scaling function obtained from the order 
parameter describes the rescaled susceptibility perfectly, without any additional adjustments of parameters. At the 
same time, this result is a remarkable validation of the functional RG approach to scaling. The critical long-range 
fluctuations, which are responsible for the critical scaling behavior, are correctly included in the effective potential 
from which we calculate the obscrvables. 

The large scaling corrections we observe for large values of the symmetry-breaking field could have implications 
for the scaling analysis of lattice QCD results, where scaling behavior is observed, but no agreement with the 0(4) 
scaling function is found. 

An obvious improvement of these results can be achieved by including a non-zero anomalous dimension in the 
calculation, work in this direction is in progress. We have also applied this approach to finite-size scaling in the 0(4) 
model a more comprehensive presentation of our results is forthcoming. 

Overall, we have obtained a result for the scaling function of the 0(4) model in three dimensions which we hope 
will prove useful, and we have demonstrated very clearly some of the remarkable relations that follow from critical 
behavior. The non-perturbative RG has proven to be a suitable tool for this application. 

Acknowledgments 

The authors would like to thank Holger Gies, Tereza Mendes, and Jan Pawlowski for fruitful discussions. This work 
was supported by the Excellence Cluster "Structure and Origin of the Universe" and by the Natural Sciences and 
Engineering Research Council of Canada (NSERC). TRIUMF receives federal funding via a contribution agreement 
through the National Research Council of Canada. 



[1] K. G. Wilson, Phys. Rev. B4, 3174 (1971). 

[2] K. G. Wilson, Phys. Rev. B4, 3184 (1971). 

[3] F. J. Wegner and A. Houghton, Phys. Rev. A8, 401 (1973). 

[4] K. G. Wilson and J. B. Kogut, Phys. Rept. 12, 75 (1974). 

[5] C. Wetterich, Phys. Lett. B301, 90 (1993). 

[6] J. Polchinski, Nucl. Phys. B231, 269 (1984). 

[7] S.-B. Liao and J. Polonyi, Phys. Rev. D51, 4474 (1995), hep-th/9403111. 

[8] S.-B. Liao, Phys. Rev. D53, 2020 (1996), hep-th/9501124. 

[9] S.-B. Liao, Phys. Rev. D56, 5008 (1997), hep-th/95 11046. 
[10] R. D. Pisarski and F. Wilczek, Phys. Rev. D29, 338 (1984). 
[11] J. Braun and H. Gies, Phys. Lett. B645, 53 (2007), hep-ph/0512085. 
[12] J. Braun and H. Gies, JHEP 06, 024 (2006), hep-ph/0602226. 

[13] Y. Aoki, G. Endrodi, Z. Fodor, S. D. Katz, and K. K. Szabo, Nature 443, 675 (2006), hep-lat/0611014. 
[14] M. D'Elia, A. Di Giacomo, and C. Pica, Phys. Rev. D72, 114510 (2005), hep-lat/0503030. 
[15] G. Cossu, M. D'Elia, A. Di Giacomo, and C. Pica (2007), arXiv:0706.4470 [hep-lat]. 
[16] D. Toussaint, Phys. Rev. D55, 362 (1997), hep-lat/9607084. 

[17] T. Schulze, J. Engels, S. Holtmann, and T. Mendes, Nucl. Phys. Proc. Suppl. 106, 498 (2002), hep-lat/01 10014. 
[18] J. Engels, S. Holtmann, T. Mendes, and T. Schulze, Phys. Lett. B514, 299 (2001), hep-lat/0105028. 
[19] F. Mendes, Braz. J. Phys. 37, 597 (2007), hep-lat/0609035. 



33 



[20] J. B. Kogut and D. K. Sinclair, Phys. Rev. D73, 074512 (2006), hep-lat/0603021. 
[21] B. Klein and J. Braun, AIP Conf. Proc. 964, 330 (2007), arXiv:0710.1038 [hep-ph]. 
[22] B. Klein and J. Braun, PoS LATTICE 2007, 198 (2007), arXiv:0710.1161 [hep-lat]. 
[23] B. Widom, J. Phys. Chem. 43, 3898 (1965). 
[24] R. B. Griffiths, Phys. Rev. 158, 176 (1967). 

[25] E. Brezin, D. J. Wallace, and K. Wilson, Phys. Rev. B7, 232 (1973). 

[26] D. J. Wallace and R. K. P. Zia, Phys. Rev. B12, 5340 (1975). 

[27] J. Engels and T. Mendes, Nucl. Phys. B572, 289 (2000), hep-lat/9911028. 

[28] A. Cucchieri and T. Mendes, J. Phys. A38, 4561 (2005), hep-lat/0406005. 

[29] R. Guida and J. Zinn-Justin, J. Phys. A31, 8103 (1998), cond-mat/9803240. 

[30] D. F. Litim and J. M. Pawlowski (1998), hep-th/9901063. 

[31] C. Bagnuls and C. Bervillier, Phys. Rept. 348, 91 (2001), hep-th/0002034. 

[32] J. Berges, N. Tetradis, and C. Wetterich, Phys. Rept. 363, 223 (2002), hep-ph/0005122. 

[33] J. Polonyi, Central Eur. J. Phys. 1, 1 (2003), hep-th/0110026. 

[34] J. M. Pawlowski (2005), hep-th/0512261. 

[35] H. Gies (2006), hep-ph/0611146. 

[36] B.-J. Schaefer and J. Wambach (2006), hep-ph/0611191. 
[37] B. Delamotte (2007), cond-mat/0702365. 

[38] J. Braun, B. Klein, and H. J. Pirner, Phys. Rev. D71, 014032 (2005), hep-ph/0408116. 

[39] J. Braun, B. Klein, and H. J. Pirner, Phys. Rev. D72, 034017 (2005), hep-ph/0504127. 

[40] J. Braun, B. Klein, H. J. Pirner, and A. H. Rezaeian, Phys. Rev. D73, 074010 (2006), hep-ph/0512274. 

[41] N. Tetradis and C. Wetterich, Nucl. Phys. B422, 541 (1994), hep-ph/9308214. 

[42] J. Braun and B. Klein (2008), (forthcoming). 

[43] B.-J. Schaefer and H.-J. Pirner, Nucl. Phys. A660, 439 (1999), nucl-th/9903003. 

[44] D. F. Litim and J. M. Pawlowski, Phys. Rev. D65, 081701 (2002), hep-th/0111191. 

[45] D. F. Litim and J. M. Pawlowski, Phys. Lett. B546, 279 (2002), hep-th/0208216. 

[46] D. F. Litim and J. M. Pawlowski, Phys. Rev. D66, 025030 (2002), hep-th/0202188. 

[47] D. F. Litim, Phys. Rev. D64, 105007 (2001), hep-th/0103195. 

[48] D. F. Litim and J. M. Pawlowski, Phys. Lett. B516, 197 (2001), hep-th/0107020. 

[49] O. Bohr, B. J. Schaefer, and J. Wambach, Int. J. Mod. Phys. A16, 3823 (2001), hep-ph/0007098. 

[50] S. A. Antonenko and A. I. Sokolov, Phys. Rev. E51, 1894 (1995), hep-th/9803264. 

[51] P. Butera and M. Comi, Phys. Rev. B56, 8212 (1997), hep-lat/9703018. 

[52] K. Kanaya and S. Kaya, Phys. Rev. D51, 2404 (1995), hep-lat/9409001. 

[53] H. G. Ballesteros, L. A. Fernandez, V. Martin-Mayor, and A. Munoz Sudupe, Phys. Lett. B387, 125 (1996), cond- 
mat/9606203. 

[54] M. Hasenbusch, J. Phys. A34, 8221 (2001), cond-mat/0010463. 

[55] G. Von Gersdorff and C. Wetterich, Phys. Rev. B64, 054513 (2001), hep-th/0008114. 

[56] C. Bervillier, A. Juttner, and D. F. Litim, Nucl. Phys. B783, 213 (2007), hep-th/0701172. 

[57] A. Pelissetto and E. Vicari, Phys. Rept. 368, 549 (2002), cond-mat/0012164. 

[58] D. F. Litim, Nucl. Phys. B631, 128 (2002), hep-th/0203006. 

[59] A. Pelissetto and E. Vicari, Nucl. Phys. B540, 639 (1999), cond-mat/9805317. 

[60] F. Karsch and E. Laermann, Phys. Rev. D50, 6954 (1994), hep-lat/9406008. 

[61] S. Aoki et al. (JLQCD), Phys. Rev. D57, 3910 (1998), hep-lat/9710048. 

[62] C. W. Bernard et al. (MILC), Phys. Rev. D61, 111502 (2000), hep-lat/9912018. 



